Introduction
Description of growth experiment.
Setup
Libraries
Flags/Config
Code
# Type of plate reader: select "Bioscreen" or "Clario"
plate_reader = "Bioscreen" ;
# normalization method:
# "initial" = normalize OD600 by subtracting the reading at t=0h
# "mean" = subtract the mean OD600 reading at each time point
norm_OD600 = "initial" ;
mytheme = theme (axis.text.x = element_text (size = 6 ),
axis.text.y = element_text (size = 6 ),
axis.title.x = element_text (size = 6 ),
axis.title.y = element_text (size = 6 ),
strip.text.x = element_text (size = 6 ),
legend.position = "bottom" ,
legend.text = element_text (size= 4 ),
aspect.ratio = 0.5 );
color_vals = c ("0" = "#BEBEBE" , "1.5" = "#68A2AD" , "2" = "#EEBC37" );
File IO
Code
# input data
file.in = "growthData_Bioscreen.xlsx" ;
path.in = "/home/tolonen/Github/actolonen/Public/Lab_Public/Growth/Microplates/Data" ;
data.in = paste (path.in, file.in, sep = "/" );
Functions
Functions to organize data
These functions take data from a specific type of plate reader and returns a data.frame with the following columns:
Hours
Well in microtiter plate
OD600
Strain (write “Blank” for lack of cells)
Dilution: dilution of cells into fresh medium in microtiter plate. For example, a 1/20 dilution = 0.05
Growth medium
Treatment (ie different carbon sources)
Function to organize Clario data
Code
organize_clario = function (data.in)
{
# output is data.frame = [Hours, Well, OD600, Treatment]
plate.map = read_excel (data.in, sheet = "Informations" , col_names = TRUE , skip = 3 );
growth.data = read_excel (data.in, sheet = "Raw data" , col_names = TRUE , skip = 0 );
# plate map variables
first.well = "A01" ;
last.well = "H12" ;
# parse plate.map
plate.map = plate.map %>%
filter (! Treatment == "Empty" );
# parse growth.data
growth.data.long = growth.data %>%
select (- Well) %>%
rename ("Time" = "...2" ) %>%
slice (- 1 ) %>% # delete first row
slice (- 1 ) %>%
mutate (Time = if_else (grepl ("min" , Time), true = Time, false = paste (Time, "0 min" , sep = " " ))) %>%
# add minutes field where lacking
separate (Time, c ("my.hours" , "hours.label" , "my.min" , "min.label" ), sep = " " ) %>%
mutate (Hours = as.numeric (my.hours) + as.numeric (my.min)/ 60 ) %>%
#mutate(Seconds = as.numeric(Seconds)) %>%
#mutate(Hours = Seconds/3600) %>%
dplyr:: select (Hours, first.well: last.well) %>%
pivot_longer (cols = first.well: last.well, names_to = "Well" , values_to = "OD600" );
growth.data.long$ OD600 = as.numeric (growth.data.long$ OD600);
# add treatment info
growth.all = inner_join (growth.data.long, plate.map, by = "Well" );
growth.all = growth.all %>%
select (Well, Hours, OD600, Strain, Medium, Treatment, Dilution)
return (growth.all)
}
Function to organize Bioscreen data
Code
organize_bioscreen = function (data.in)
{
# output is data.frame = [Hours, Well, OD600, Treatment]
plate.map = read_excel (data.in, sheet = "Informations" , col_names = TRUE , skip = 29 );
growth.data = read_excel (data.in, sheet = "Raw data" , col_names = TRUE , skip = 2 );
# plate map variables
first.well = "1" ;
last.well = "199" ;
# parse plate.map
plate.map = plate.map %>%
filter (! Treatment == "Empty" ) %>%
mutate (Well = as.character (Well));
# calculate elapsed time in hours
time1 = hms ("00:00:00" ); # start time
growth.data = growth.data %>%
mutate (Time_hms = hms (Time)) %>%
mutate (Time_Elapsed = Time_hms - time1) %>%
mutate (Hours = hour (Time_Elapsed) +
(minute (Time_Elapsed) / 60 ) +
(second (Time_Elapsed)/ 3600 ));
growth.data.long = growth.data %>%
select (- Blank) %>%
dplyr:: select (Hours, first.well: last.well) %>%
pivot_longer (cols = first.well: last.well, names_to = "Well" , values_to = "OD600" );
growth.data.long$ OD600 = as.numeric (growth.data.long$ OD600);
# remove contaminated wells
#growth.data.long = growth.data.long %>%
# filter(!Well %in% c("A04", "A05", "A06", "C03", "E02", "G01", "G04", "G05"));
# add treatment info
growth.all = inner_join (growth.data.long, plate.map, by = "Well" );
growth.all = select (growth.all, Well, Hours, OD600, Strain, Medium, Treatment, Dilution);
return (growth.all)
}
Functions to normalize OD600 data
Function to normalize OD600 readings by subtracting OD600 at time zero
Code
# Function to normalize OD600 readings by subtracting the OD600 at t=0h
norm_initial = function (growth.cells)
{
# find first timepoint
min_time = min (growth.cells$ Hours);
# calc initial OD600
initial.OD = growth.cells %>%
filter (Hours == min_time) %>%
rename (OD600_t0 = OD600) %>%
select (Well, OD600_t0);
# normalize OD600 by subtracting initial reading
growth.cells = left_join (growth.cells, initial.OD, by = "Well" );
growth.cells = growth.cells %>%
mutate (OD600 = OD600 - OD600_t0) %>%
select (- OD600_t0);
return (growth.cells);
}
Function to normalize OD600 readings by subtracting mean OD600 at the same time point
Code
# substract blank mean from OD measurements
norm_mean = function (growth.blanks, growth.cells)
{
# calc mean blank
growth.blanks = growth.blanks %>%
group_by (Hours)%>%
mutate (OD600_blank = mean (OD600)) %>%
ungroup;
growth.blanks = growth.blanks %>%
select (Hours, OD600_blank) %>%
distinct ();
# subtract mean blank
growth.cells = left_join (growth.cells, growth.blanks, by = c ("Hours" ));
growth.cells = growth.cells %>%
mutate (OD600_norm = OD600 - OD600_blank);
return (growth.cells);
}
Function to run growthcurver
Code
# function to run Growthcurver analysis.
# input = data.frame(Hours, Treatment, OD600_mean, OD600_sd)
# output =
run_growthcurver = function (growth.treatment)
{
# select cols of interest
growth.fit = growth.treatment %>%
select (Hours, OD600_mean) %>%
rename (time = Hours) %>% # growthcurver requirement
mutate (blank = 0 );
# calc growth rate information
gc_fit = SummarizeGrowthByPlate (growth.fit, bg_correct = "blank" );
gc_fit = gc_fit %>%
slice (- 2 ); # get rid of blank data
# components of SummarizeGrowthByPlate object
k_maxOD = gc_fit$ k; # carrying capacity
n0_minOD = gc_fit$ n0; # initial OD600
r_intrinsicGrowth = gc_fit$ r; # intrinsic growth rate
genTime_hours = gc_fit$ t_gen; # generation time
auc_areaCurve = gc_fit$ auc_l; # AUC, integral of logistic eq
auc_areaTrapezoid = gc_fit$ auc_e; # AUC, area of trapezoid
sigma_fit = gc_fit$ sigma; # goodnesss of fit
my_note = gc_fit$ note; # note if poor fit
# make data.frame of fit data
Treatment = growth.treatment$ Treatment[1 ];
fitdata = data.frame (Treatment,
k_maxOD,
n0_minOD,
r_intrinsicGrowth,
genTime_hours,
sigma_fit,
my_note);
return (fitdata);
}
Main
Organize data
Code
if (plate_reader == "Clario" )
{
growth.all = organize_clario (data.in);
}
if (plate_reader == "Bioscreen" )
{
growth.all = organize_bioscreen (data.in);
}
growth.short = head (growth.all, 5 );
table.growth = kable (growth.short, caption = "Table: tidy growth data" );
table.growth
Table: tidy growth data
1
0.0347222
0.190
2
GS2
0
0,05
2
0.0347222
0.198
2
GS2
0
0,05
3
0.0347222
0.195
2
GS2
0
0,05
4
0.0347222
0.188
2
GS2
1
0,05
5
0.0347222
0.190
2
GS2
1
0,05
Plots
Plot blanks
Code
growth.blanks = growth.all %>%
filter (grepl ("Blank" , Strain));
plot.blanks = ggplot (growth.blanks, aes (x= Hours, y= OD600, group= Well)) +
geom_point (size= 1 , color = "green" ) +
geom_line (linewidth= 0.1 , color= 'black' ) +
theme_bw () +
xlab ("Hours" ) +
ylab ("OD(600)" ) +
coord_cartesian (
xlim = c (0 , max (growth.blanks$ Hours)),
ylim = c (0 , max (growth.blanks$ OD600) + 0.2 ))+
facet_wrap (~ Treatment, ncol= 3 )+
scale_x_continuous (breaks= seq (0 , max (growth.blanks$ Hours), 12 ))+
scale_y_continuous (breaks= seq (0 , max (growth.blanks$ OD600 + 0.1 ), 0.1 ))+
theme_classic ()+
mytheme;
grid.arrange (plot.blanks, ncol= 1 );
Normalize OD600 data: subtract blanks from wells with cells.
Code
# focus on well with cells
growth.cells = growth.all %>%
filter (! grepl ("Blank|Empty" , Strain));
if (norm_OD600 == "mean" )
{
growth.cells = norm_mean (growth.blanks, growth.cells);
}
if (norm_OD600 == "initial" )
{
growth.cells = norm_initial (growth.cells);
}
Plot growth: for each strain (panel) compare growth in each treatment (lines)
Code
plot.cells = ggplot (growth.cells, aes (x= Hours, y= OD600, group= Well, color= Treatment)) +
geom_point (size= .05 ) +
geom_line (linewidth= 0.1 , color= 'black' ) +
theme_bw () +
xlab ("Hours" ) +
ylab ("OD(600)" ) +
coord_cartesian (
xlim = c (0 , max (growth.cells$ Hours)),
ylim = c (- 0.05 , max (growth.cells$ OD600) + 0.1 ))+
facet_wrap (~ Strain, ncol= 7 )+
scale_x_continuous (breaks= seq (0 , max (growth.cells$ Hours), 24 ))+
scale_y_continuous (breaks= seq (0 , max (growth.cells$ OD600 + 0.1 ), 0.4 ))+
scale_color_manual (values = c (color_vals))+
theme_classic ()+
mytheme;
grid.arrange (plot.cells, ncol= 1 );
Calculate treatment means for each strain
Code
growth.cells = growth.cells %>%
group_by (Hours, Treatment, Strain) %>%
mutate (OD600_mean = mean (OD600)) %>%
mutate (OD600_sd = sd (OD600)) %>%
ungroup ();
Plot growth: for each treatment (panel) compare growth of each strain (lines)
Code
plot.cells.means = ggplot (growth.cells, aes (
x= Hours,
y= OD600_mean,
group = Strain,
ymin = OD600_mean - OD600_sd,
ymax = OD600_mean + OD600_sd)) +
geom_point (size= .05 , color = "blue" ) +
geom_line (linewidth= 0.1 , color= 'black' ) +
# geom_errorbar(width = 0, color = "blue")+
theme_bw () +
xlab ("Hours" ) +
ylab ("OD(600)" ) +
coord_cartesian (
xlim = c (0 , max (growth.cells$ Hours)),
ylim = c (- 0.05 , max (growth.cells$ OD600) + 0.1 ))+
facet_wrap (~ Treatment, ncol= 7 )+
scale_x_continuous (breaks= seq (0 , max (growth.cells$ Hours), 24 ))+
scale_y_continuous (breaks= seq (0 , max (growth.cells$ OD600 + 0.1 ), 0.4 ))+
theme_classic ()+
mytheme;
grid.arrange (plot.cells.means, ncol= 1 );
Plot growth means: for each strain (panel) compare growth in each treatment (lines)
Code
plot.cells_strains = ggplot (growth.cells, aes (
x= Hours,
y= OD600_mean,
group= Treatment,
color= Treatment))+
#ymin = OD600_mean - OD600_sd,
#ymax = OD600_mean + OD600_sd)) +
geom_point (size= .05 ) +
geom_line (linewidth= 0.1 , color= 'black' ) +
#geom_errorbar(width = 0, color = "blue")+
facet_wrap (~ Strain, ncol= 7 )+
theme_bw () +
xlab ("Hours" ) +
ylab ("OD(600)" ) +
coord_cartesian (
xlim = c (0 , max (growth.cells$ Hours)),
ylim = c (- 0.05 , max (growth.cells$ OD600) + 0.1 ))+
scale_x_continuous (breaks= seq (0 , max (growth.cells$ Hours), 12 ))+
scale_y_continuous (breaks= seq (0 , max (growth.cells$ OD600 + 0.1 ), 0.4 ))+
theme_classic ()+
mytheme;
grid.arrange (plot.cells_strains, ncol= 1 );
Calc max OD for each strain in each treatment
Code
# focus on data of interest
growth_max = growth.cells %>%
dplyr:: select (OD600_mean, Strain, Treatment) %>%
distinct ();
# calc max OD for each strain in each treatment
growth_max = growth_max %>%
group_by (Strain, Treatment) %>%
mutate (OD600_max = max (OD600_mean)) %>%
ungroup %>%
dplyr:: select (Strain, Treatment, OD600_max) %>%
distinct ();
Plot max OD: for each strain in each treatment
Code
# remove WT
growth_max = growth_max %>%
filter (! Strain == "WT" );
plot_maxOD = ggplot (growth_max, aes (
x = Treatment,
y = OD600_max,
color = Treatment))+
geom_boxplot (outlier.size = 0 )+
geom_jitter (aes (text = Strain), position= position_jitter (0.2 ))+
xlab ("Treatment" )+
ylab ("max(OD600)" )+
coord_cartesian (ylim= c (0 , 0.6 ))+
scale_y_continuous (breaks= seq (0 , 0.6 , 0.1 ))+
scale_color_manual (values = c (color_vals))+
theme_classic ()+
mytheme;
plotly_maxOD = ggplotly (plot_maxOD, tooltip= "text" );
plotly_maxOD
Fig 5: max OD600 for each strain in each treatment.
Quantify growth rates using growthcurver
Select treatments for growthcurver analysis
Code
growth_gc = growth.cells %>%
filter (Strain %in% c ("2" )) %>%
filter (Treatment %in% c ("0" , "1" ));
Select log phase growth data for growthcurver
Code
# data.frame of maxODs for each treatment
data.maxOD = growth_gc %>%
group_by (Treatment) %>%
mutate (maxOD = max (OD600_mean)) %>%
filter (OD600_mean == maxOD) %>%
ungroup;
# filter for first time OD reaches max OD
data.maxOD = data.maxOD %>%
group_by (Treatment) %>%
arrange (Hours) %>%
filter (row_number ()== 1 ) %>%
ungroup () %>%
mutate (time_maxOD = Hours) %>%
select (Treatment, time_maxOD);
growth.cells.means.log = left_join (growth_gc, data.maxOD, by = "Treatment" );
# filter for data in log phase
growth.cells.means.log = growth.cells.means.log %>%
filter (Hours < (time_maxOD + 3 )) %>%
select (Hours, Treatment, OD600_mean);
Run growthcuver
Code
# get list of treatments/strains
treatments = unique (growth.cells.means.log$ Treatment);
# declare empty data.frame for fit data
fit.data.all = data.frame (Treatment = character (),
k_maxOD = double (),
n0_minOD = double (),
r_intrinsicGrowth = double (),
genTime_hours = double (),
sigma_fit = double (),
my_note = character ());
# run growthcurver on each treatment/strain
for (my.treatment in treatments)
{
data.fit = growth.cells.means.log %>%
filter (Treatment == my.treatment);
fit.data = run_growthcurver (data.fit);
fit.data.all = rbind (fit.data.all, fit.data);
}
# add fit data to growthdata
growth.cells.fit = left_join (growth.cells.means.log,
fit.data.all,
by= "Treatment" );
growth.cells.fit = growth.cells.fit %>%
mutate (my_fit = k_maxOD / (1 + (((k_maxOD - n0_minOD) / n0_minOD) * exp (1 )^- (r_intrinsicGrowth * Hours))));
Plot growthcurver data
Code
# subsample data for plotting
growth.cells.fit.plot = growth.cells.fit %>%
filter (round (Hours, 2 ) %% 2 == 0 );
fitplots = ggplot (growth.cells.fit.plot, aes (x= Hours,
y= OD600_mean))+
facet_wrap (~ Treatment)+
#geom_line(size=0.2, color='black') +
geom_point (size= 1 ) +
geom_line (aes (x= Hours, y= my_fit, color = 'red' ))+
xlab ("time" ) + ylab ("OD(600)" )+
coord_cartesian (
xlim = c (0 , 40 ),
ylim = c (0 , 0.6 ))+
theme_classic ()+
mytheme;
fitplots
Table of growthcurver data
Code
growth.cells.fit.table = growth.cells.fit %>%
select (Treatment, k_maxOD, r_intrinsicGrowth, genTime_hours, sigma_fit) %>%
mutate (k_maxOD = round (k_maxOD, 2 )) %>%
mutate (genTime_hours = round (genTime_hours, 2 )) %>%
mutate (r_intrinsicGrowth = round (r_intrinsicGrowth, 2 )) %>%
mutate (sigma_fit = round (sigma_fit, 2 )) %>%
distinct ();
table1 = kable (growth.cells.fit.table, caption = "Growth curve data for each treatment. " );
table1
Growth curve data for each treatment.
0
0.45
0.87
0.80
0.01
1
0.27
0.76
0.92
0.01